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I, Task Synopsis 


Air-breathing hypersonic vehicles, with fully integrated airframe and propulsion systems, 
are expected to face the most complicated design, trajectory, and control coupling challenges 
of any class of flight vehicles ever considered. Assessment of appropriate trajectory management 
and control strategies must be accomplished during the conceptual phases of vehicle design to 
assure a high probability that the fully integrated system will accomplish mission objectives. Also, 
a valid mathematical basis is needed to rate the relative performance of alternative vehicle design 
concepts and/or subsystem technologies against a common mission goal. To address these 
concerns, a generically applicable methodology has been developed to construct trajectories 
and control strategies while accommodating physically derived constraints on vehicle operations. 

This task has adapted a Draper Laboratory developed tool which generates optimal 
trajectory/control histories in an integrated manner to the treatment of single-stage-to-orbit 
air-breathing hypersonic vehicles. The mathematical basis for implementing the methodology 
as a two-point boundary value problem solution technique is documented. The tool formulation 
permits an assessment of the entire near-minimum-fuel trajectory and desired control strategy 
from takeoff to orbit while satisfying physically derived inequality constraints and while achieving 
efficient propulsive mode phasing. A strategy has also been formulated and is documented to 
construct preliminary trajectory and control history representations with less computational 
burden than required for the overall flight profile assessment (by a partitioning of the trajectory 
into boundary condition matched segments). 

A tabulated data version of an example hypersonic vehicle model has been used in a 
demonstration of the integrated analysis methodology. To assure good numerical behavior of 
the algorithm (no discontinuities in the interpolated values or their first derivatives), it has been 
necessary to develop a smoothing routine for the tabulated data. A cubic spline curve smoothing 
routine that can function with any multidimensional tabulated vehicle model was developed 
implemented, and documented. ’ 

A demonstration of the analysis methodology was accomplished using the example vehicle 
model and the data smoothing routine. A constrained near-fuel-optimal trajectory was developed 
from horizontal takeoff to 20,000 ft/sec relative air speed while aiming for a polar orbit. The 
inequality constraint imposed on the system was a dynamic pressure limit of 1 000 psf. 

The demonstration case provided examples of the kind of information that the analysis 
methodology can develop for a wide range of HSV designs and mission applications. A 
near-fuel-optimal trajectory that satisfies the applied constraints was determined. Previously 
unspecified propulsive discontinuities were located. Flight regimes demanding rapid attitude 
changes were identified, dictating control effector and closed-loop G&C design requirements 
based on the required control response. Available closed-loop controller authority was 
ascertained after evaluating effector use for vehicle trim. Also, inadequacies in vehicle model 
representations and specific subsystem models with insufficient fidelity were determined (as 
manifested by unusual control characteristics and/or excessive sensitivity to uncertainty). 

After completion of analysis of the information resulting from use of the integrated 
trajectory/control methodology, closed-loop guidance and control design requirements can be 
established. This information must precede specification of on-board advanced information 
processing system design characteristics. 
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III. List of Symbols and Acronyms 


a 

AIPS 

B 


c 

C 

C n 



Cj 


C, 

Cl. 

Cl. 


c 

c 


rn 


C 


c 


Cq 

c r 


the atmospheric speed of sound (a function of r) 
advanced information processing system 
a matrix term in the trajectory optimization variation equations 
the vehicle longitudinal reference length 
a vector term in the variation equations 
the vehicle drag coefficient 

the drag increment coefficient for the basic vehicle derived from the Langley data 
base (a function of « and M ,) 

the drag increment coefficient of the combined left and right elevons derived from 
the Langley data base (a function of a , 6 P , and M ,) - it has a value of zero after the 
rocket engine turn-on switch time 

a coefficient to establish the desired weight given to cost improvement contributions 
to variation terms 

a coefficient to establish the desired weight given to cost improvement contributions 

of variation terms on the i^ trajectory leg 
the vehicle lift coefficient 

the lift increment coefficient for the basic vehicle derived from the Langley data base 
(a function of « and Af ,) 

the lift increment coefficient of the combined left and right elevons derived from the 
Langley data base (a function of a , 6 P , and M f ) - it has a value of zero after the 
rocket engine turn-on switch time 

the vehicle pitch moment coefficient relative to the moment reference center 
the pitch moment increment coefficient for the basic vehicle, relative to the moment 
reference center, derived from the Langley data base (a function of ot and M f ) 
the pitch moment increment coefficient of the combined left and right elevons, rela- 
tive to the moment reference center, derived from the Langley data base (a function 
of a , 6 e , and M y ) - it has a value of zero after the rocket engine turn-on switch time 

the vehicle pitching moment increment coefficient with pitch rate derivative, relative 
to the moment reference center, derived from the Langley data base (a function of a 
and M ,) 

the dynamic pressure constraint cost weighting term 

the vehicle air-breathing engine thrust coefficient derived from the Langley data 
base with dimension ft 2 (a function of 4> a , M , , and weakly of Q) - it has a value of 
zero after the air-breathing engine turn-off switch time 

a diagonal matrix of constant coefficients to establish the desired weights given to 
equality constraint violation improvements in variation terms 
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C y a diagonal matrix of constant coefficients to establish the desired weights given to 

equality constraint violation improvements in variation terms on the i^ trajectory 
segment 

D the drag force acting on the vehicle 

DFY Draper fiscal year 

/ the derivative of the state vector .v with respect to time 

/ , the derivative of x , with respect to time 

f R The contribution to the equation for f due to the rocket thrust 

f s ( the derivative of the state vector x with respect to time evaluated at the jth 

switch point 

f u the partial derivative of the vector / with respect to the control vector u 

f x the partial derivative of the vector f with respect to the state vector x 

F the LVLH vector of forces acting on the vehicle 

F g the gravity force acting on the vehicle 

F T the force component parallel to the vehicle free stream velocity vector 

F N the force component normal to the vehicle free stream velocity vector, pointing 

down in the vehicle plane of symmetry 

F i the force component normal to F 2 and F 3 , pointing toward increasing p 

F 2 the force component parallel to the equatorial plane, pointing eastward, 

normal to F 3 

F 3 the force component parallel to r - points toward Earth 
g a vector term in the variation equations 

Sfo a constant representing the number of pounds mass per slug (approx. 32.2) 

G&C guidance and control 

h vehicle altitude 

II the Hamiltonian 

H u the partial derivative of the Hamiltonian H with respect to the control vector u 
HSV hypersonic vehicle 

I j j a scalar influence function of cost 

* j.j, a scalar influence function of cost for the i^ 1 trajectory leg 

I SPa the air-breathing engine specific impulse derived from the Langley data base with 
dimension sec (a function of 4> a , M f , and weakly of Q) 

1 S p K the rocket engine specific impulse (a constant) 

I yy the vehicle pitch moment of inertia derived from the Langley data base with 
dimension slug-ft 2 

I y j a vector influence function of constraints and cost 

I y _ J( a vector influence function of constraints and cost for the i th trajectory segment 

I y y a matrix influence function of constraints 

I y ( y 4 3 matrix influence function of constraints for the i^ trajectory segment 
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IR&D 

J 

J, 

Js 

L 

L u 

L x 

K Ql 

K. 


K, 
K 3 

L 

L 

L, 


LaRC 

LVLH 

m 
m a 
m R 
m o 
M 
M, 
M q 

NM 

psf 

q 

Qo 

Q 

Q D 

r 

r e 


internal research and development 
a mathematical cost imposed on the system 

a mathematical cost imposed on the system for the i th ascent trajectory leg 
a heuristically determined specified cost improvement 
the distributed mathematical cost term 

the partial derivative of the distributed mathematical cost term L with respect to the 
control vector u 

the partial derivative of the distributed mathematical cost term L with respect to the 
state vector x 

the distributed cost weighting term on the i^ 1 ascent trajectory leg 
a coefficient in the equation used to find 6 e for vehicle moment trim - a function of a 
and m 

a coefficient in the equation used to find 6 e for vehicle moment trim - a function of a 
and m 

a coefficient in the equation used to find 6 e for vehicle moment trim - a function of 
a , m , and M f 

the distributed mathematical cost term 
the lift force acting on the vehicle 

the distributed mathematical cost term for the i^ 1 ascent trajectory leg 

the partial derivative of the distributed mathematical cost term L, with respect to the 

control vector x 

the partial derivative of the distributed mathematical cost term L with respect to the 

control vector u 

Langley Research Center 

the local vertical/local horizontal reference frame 

the vehicle mass 

the air-breathing engine fuel mass 
the rocket engine propellant mass 
the vehicle mass at takeoff 

a matrix term used in the trajectory optimization variation equations 

the free stream Mach number 

the pitch moment acting on the vehicle 

nautical miles 
pounds per square foot 
the vehicle pitch rate 
the vehicle takeoff pitch rate 
the dynamic pressure 

the desired dynamic pressure constraint bound 
the vehicle distance from the center of the Earth 
the radius of the Earth 
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r larg the target circular orbit altitude 

r o the vehicle distance from the center of the Earth at takeoff 

S ref the vehicle reference wing area 

SSTO single-stage-to-orbit 

t time 

l s the switch time vector 

t s the jth switch time 

i 

T the thrust force acting on the vehicle 

T R the rocket thrust at full throttle - it has a value of zero before the rocket engine 
turn-on switch time 
a the control vector 

u a the control for the air-breathing engine 

u, the ifr* element of the control vector 

u R the control for the rocket 

u 0 a step function which changes from 0 to 1 at a functional value designated in 

brackets 

U a diagonal matrix of functions of time to weight different elements of the control 

variation vector 

v f the vehicle velocity magnitude relative to the free stream 

Vtarg the target velocity magnitude at the desired circular orbit condition 

Utarg, the target velocity magnitude for the end of i^ ascent trajectory leg 

v i the free stream velocity component normal to v 2 and u 3 , pointing toward 

increasing p 

v , o the northerly component of the takeoff velocity 

v 2 the free stream velocity component parallel to the equatorial plane, pointing 

eastward, normal to v 3 

v 2q the easterly component of the takeoff velocity 

v 3 the vehicle free stream velocity component parallel to r (points toward Earth) 

v 3q the vehicle free stream velocity component at takeoff parallel to r (points toward 

Earth) 

V a diagonal matrix of constant coefficients to weight different elements of the switch 

time variation vector 
x the state vector 

x cc the longitudinal distance from the current vehicle center of gravity to the moment 
reference center derived from the Langley data base with dimension ft 
x, the ith component of the state vector x 

x R the longitudinal distance from the vehicle moment reference center to the point of 

rocket thrust application 
x o the vehicle state vector at takeoff 

a the vehicle angle of attack 
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A wide variety of HSV design concepts and applicable subsystem technologies are under 
consideration. Unique trajectory management strategies are likely to be necessary for each 
proposed vehicle configuration to assure efficient payload performance while satisfying critical 
physical constraints. The generic analysis technique implemented under this subtask permits a 
common tool to be applied to all proposed HSV SSTO configurations. Use of the tool supplies 
information needed for a valid comparison of expected performance and key sensitivities for 
alternative vehicle concepts. 

The methodology applies a generalized gradient two-point boundary value problem solution 
algorithm to the HSV problem. Equality constraints are adjoined to the problem. Inequality 
constraints are treated by either mathematical remapping of the function space or penalty 
functions. Dynamic discontinuities, such as phasing of propulsive modes, are treated by an 
extension of the control vector to include the times at which the discontinuities occur. The 
integrated analysis data are generated by optimizing a mathematical performance index which 
emphasizes payload performance. 

iv.B. Longer Term Research Objective s 

The task statement of work identifies several investigation paths which are intended as 
directions for research under possible continuations of the effort. These are summarized below: 

IV.B.1. Treatment of Vehicle Thrust Direction and Pitch Moment Changes with Angle of 
Attack 

The initial example HSV model included a propulsion system that represented thrust 
magnitude as independent of angle of attack and forced thrust direction to always be in the 
lonqitudinal body axis. More realistic vehicle models specify the thrust magnitude and direction 
as functions of angle of attack, which also implies associated longitudinal moment variations. 
This subtask would characterize the changes in HSV trajectory management and control 
strategies that result from explicit consideration of angle of attack effects on propulsion 
performance. 

IV.B.2. Incorporation of Thermal Constraints Into Trajectory Results 

The integrated trajectory/control analysis methodology supports treatment of inequality 
constraints which are functionally dependent on both states and controls. However, the initial 
demonstration was restricted to studying the effects of dynamic pressure inequality constraints 
which are only state dependent. Treatment of thermal constraints would address an important 
physical performance limitation that has both state and control dependencies. This subtask would 
incorporate a vehicle stagnation point thermal flux constraint and assess its qualitative effects on 
mission performance and G&C strategies. 


IV.B.3. HSV G&C Performance Robustness Studies 

Uncertainty in the environment model and expected vehicle dynamics are important design 
considerations for HSV G&C systems. The atmosphere can experience significant density 
variations at the higher HSV air-breathing altitudes. The aerodynamic and propulsion performance 
of HSVs may not be accurately modelled before there is significant flight experience and may 
vary between vehicles or for a given vehicle on different flights. The likelihood of successfully 
accomplishing specific mission objectives depends on understanding the effects of these 
uncertainties and designing a G&C system that is robust enough to handle them. This subtask 
would evaluate the performance of the near-fuel-optimal HSV trajectories and candidate 
closed-loop G&C systems in the presence of aerodynamic variations, propulsion modelling 
uncertainties, and atmospheric disturbances. These uncertainty effects would be evaluated for 
a range of applicable inequality constraint bounds. 
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V, Progress for June 1990 to March 1991 
V.A, Formulation of an Integrated Traiectory/Control Analysis Methodology 

V.A.1. Overview 

A major objective of this task has been to demonstrate the ability of the Draper integrated 
trajectory/control analysis methodologies to derive useful G&C design information when applied 
to realistic HSV ascent models. To accomplish this, the resulting tools were adapted to develop 
near-minimum-fuel, horizontal take-off-to-orbit trajectories for an example configuration a 
winged-cone vehicle which is defined in atabulated parametric database [1 1 ]. The model features 
lift, drag, side forces, and moments as functions of a , Mach number, and 6 e . Also, the air-breathing 

propulsion system thrust and specific impulse are represented as a function of a , Mach number, 
and Q . 

The target orbit of interest in the study is polar, circular, and low altitude (e.g. 1 1 0 NM). The 
target requires some rocket thrust outside the sensible atmosphere to achieve the circularization. 
Therefore, full trajectory development requires evaluation of both air-breathing engine and rocket 
throttle setting histories. This requires treatment of the flight condition at which the air-breathing 
engine is fully shut down, the flight condition at which the rocket engine is first turned on, and 
aerosurface and rocket engine deflection histories. A suitable trajectory must minimize fuel 
consumption, satisfy physical design limits such as an upper bound on dynamic pressure and 
must match both desired takeoff conditions and the desired target orbit state. 

An open loop trajectory optimization algorithm general enough to handle all the 
considerations of concern for HSVs has already been developed and demonstrated for other 
HSV designs [2,4,1 2]. While the mathematics of the algorithm remain similar in each application, 
it is necessary to develop a performance index and apply constraints unique to each application’ 
In addition, each vehicle model and associated partial derivatives must be defined. These features 
as well as a standard atmosphere model must be implemented in a form compatible with the 
optimization algorithm. 

Air-breathing HSVs have highly nonlinear dynamics during powered flight that vary greatly 
with Mach number. In addition to extensive coupling between propulsion and control, there is 
little relevant flight experience. Construction of nearly fuel-optimal take-off-to-orbit trajectories 
cannot be accomplished from existing analytic formulas or by reliance on previous results 
However, use of an applicable optimization algorithm can generate a desired trajectory. The 
algorithm must explicitly consider an appropriate performance index, boundary conditions 
physically derived inequality constraints, and the phasing of discrete propulsion system 
operations. 7 

The following specific methodology capabilities have been demonstrated in the trajectory 
development for the example HSV: 

- Use of the example vehicle model's multidimensional tabulated database. 

- Horizontal launch to near orbital velocities on a polar orbit track. 

- Application of an upper bound on dynamic pressure. 

- Achievement of near minimum propellant consumption. 

- Consideration of elevon aerodynamics and longitudinal rotation dynamics. 

- Development of an algorithm feature that can determine nearly optimal phasing of the 
air-breathing engine shutdown and rocket ignition. 
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V.A.2. Data Smoothing Techniques 

For the purpose of interpolating the aerodynamic and propulsion data, a new multidi- 
mensional cubic spline algorithm was developed to replace the B-spline method that has been 
used previously at Draper for similar applications. 

B-splines are a basis set of overlapping polynomial splines from which interpolated values 
are obtained through linear combination over each interval. The main advantage of this approach 
is its generality. B-splines can be applied to multidimensional interpolation problems and can be 
used to obtain any desired degree of smoothness (continuous N" 1 order derivatives) at any point 
(including extreme points) in the data set. Continuity of the derivatives is achieved through the 
use of knot sequences along the independent axes. For more details about B-spline algorithms 
and their properties refer to ref. [13]. The disadvantage of using B-splines is that for problems 
not requiring such flexibility there is a significant degree of overhead in terms of transforming the 
data into B-spline tables, performing tensor algebra for multidimensional problems, and 
generating the required software. 

The new interpolation algorithm is a multidimensional extension of the basic cubic spline 
procedure found in many calculus textbooks. It is also a variation of the so called "Cubic Hermite 
Interpolation" and provides a very efficient local scheme that is continuous up to the second-order 
derivatives. For problems requiring only smooth first-order derivatives it is a much faster algorithm 
than the B-spline approach. 

The data set is assumed to be in the form of an N-dimensional table with uniformly spaced 
grid points along each of the axes of the independent variables. Prior to performing any inter- 
polations, the entire data set is preprocessed to generate derivative information at the grid points. 
Since no a priori information is available to bias the choice of these derivatives, they are computed 
based on the straight line slope between the bounding grid points on each axis. For example, 
in one dimension, along the x-axis, x being an independent variable, at grid point *( 2 ): 


df _ /(-v(3))~ nx(JJ) 
rtx a ( 3 ) - a* ( 1 ) 


(Eq. V.1) 


At the data set boundaries the last interval slope is used. For example, at x(0): 


df _ /•(-v(l))-/(-v(0)) 

dx .v ( I ) - .v ( 0 ) 


(Eq. V.2) 


Precomputing the derivatives is not necessary, however it greatly improves the efficiency 
of the on-line interpolation. 

The interpolation algorithm is best visualized using a two-dimensional example as illustrated 
in figure V-1. Consider the case of a dependent variable, such as drag coefficient, tabulated as 
a function of two independent variables, such as Mach Number and angle-of-attack. During an 
interpolation, one first finds the two-dimension grid cell containing the current flight condition. In 
the figure (xcoord, ycoord) represents the flight condition (Mach number and angle of attack), 
and the coordinates (a,c), (b,c), (b,d), and (a,d) represent the bounding cell. The function values, 
and the associated derivatives (along each axis) , are then looked up for each of these coordinates. 
Using these values, several one-dimensional cubic splines are constructed in the same direction 
along the vertices of the bounding region. Note from the figure that two splines were required 
in the x-axis direction for a two-dimensional problem (labeled as fitl and fit2). In three dimensions, 
this first step would have required four splines, and in more dimensions 2™-< splines. These 
one-dimensional splines are then used to obtain function values at the coordinates (xcoord, c) 
and (xcoord.d), while the derivatives at these points are obtained by linearly interpolating the 
slopes at the end-points of each spline. This completes the first interpolation step which effectively 
reduces the problem by one dimension. Now a simple cubic spline (fit3) can be generated along 
the y-axis using the new function values and slopes. The final interpolated value at (xcoord , ycoord) 
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is then obtained from this last spline at y=ycoord. This method works equally well for hiqher 
order interpolations, using one-dimensional splines to reduce the number of dimensions 
step-by-step until one final spline gives the desired interpolated value. It can be shown that the 
mterpolated values using this algorithm are continuous with respect to the independent variables 
and that their first-order derivatives are smooth. ’ 


The interpolation of aerodynamic and propulsion data was performed using the cubic spline 
algorithm described above. Although a general N-dimensional routine could have been 
implemented, separate two-dimensional and three-dimensional versions were written for speed 
and simplicity of the code. The interpolation software is divided into three sets of routines 
representing a preprocessor to generate derivatives, a table look up function to locate the cell 
boundaries, and then the actual code to perform the interpolation. The resulting interpolation 
f° de •seff'C'ent in terms of code size as well as execution speed. The only penalty is related to 
the additional storage space required to store the derivatives. 


V.A.3. Mathematical Basis of Analysis Algorithm 

V.A.3.a. The Optimization Problem Formulation and Features 

The trajectory optimization algorithm uses a generalized steepest descent aradient 
two-point boundary value problem solution technique. The most general form of the Draper 
developed methodology is capable of solving simultaneously for vehicle configuration traiectorv 
dynamic discontinuity times, and controls [2,4], It can accommodate equality and inequality 
constraints while optimizing with respect to a performance index of states, vehicle desian 
parameters, controls, and time. In the current application, vehicle geometry is assumed fixed, 
eliminating consideration of design parameters. The optimization performance index J as applied 
to this problem has the form: 


J(T)= ())(.V(T)) 



L(-v(f ), u(t))dt 


(Eq. V.3) 


where x is the time at which the target orbit state is reached, x is the vehicle state vector: 



(Eq. V.4) 
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and u is the control vector: 


A,\ 


/«,A 


it 


u 
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(Eq. V.5) 
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The relationship between the control vector elements and the actual propulsion system 
throttle variables are: 

2 (Eq- V.6) 


<*>„ = (“«)' 

I 




I + (<*/?)' 


(Eq. V.7) 


The terminal cost term <t> is used to minimize the consumed fuel mass. 

cj> = m 0 - m(x) 


(Eq. V.8) 


The integral function I . treats inequality constraints. In this application where dynamic 
pressure is constrained, I is a function only of states since dynamic pressure is a function of 
velocity states and density, and density is a function of r when using a standard atmosphere 
model. 

l.-CVQ-Q^uJQ-Qol (Eq - V9) 


where: 


<?- 


pO/) : 

2 


(Eq. V.10) 


and where it 0 [ Q - Q D j is a unit step function that changes from zero to one at Q - Q D . 

A costate like function A is adjoined to the problem to explicitly treat a vector of applicable 
equality constraints functions called V. The specific equations for A are developed in section 
V B 2 c in Eqs. V.45 and V.46. With the initial take-off conditions given, the equality constraints 
are needed to enforce the target orbit conditions at time x . One smoothly varying target orbit 
state equality constraint is needed to determine x . The target orbital velocity is used in this 
application and is therefore excluded from the adjoined function list. This trajectory integration 
cutoff function called Q helps determine boundary conditions on the costate function X. , while 
the remaining equality constraints determine boundary conditions used in the adjoined function 
\ The specific formulation for X is also detailed in section V.B.2.C in Eqs. V.42,43. 

An algorithm feature that locates the most efficient times for dynamic discontinuity events, 
l is used to optimize propulsive phasing of the air-breathing engine and rocket. The 

discontinuities are factored into the vehicle state dynamics model (see Eq. V.38) and are treated 
^ the cosSteant? adjoined function dynamics (see Eqs. V.47.48). A discrete a,r-breath,ng engine 
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shutdown time and a separate discrete rocket ignition time can be evaluated by this means. The 
implementation permits staggered or overlapping engine operation as dictated by minimum fuel 
considerations. 

Given an initial control history and propulsive phasing time guess, the algorithm iteratively 
computes perturbations to the controls and discontinuity times to improve the performance index 
while maintaining boundary conditions. Additive perturbation contributions are derived to improve 
the cost as well as to drive down violations of the equality constraints. Weighting functions balance 
the relative contributions of each additive perturbation term. 

Numerical stability of the algorithm requires that the initial trajectory guess nearly satisfies 
the equality constraints even if far from optimal in performance. Since the vehicle dynamics are 
very dissimilar from a simple rocket, a strategy to develop a suitable initial guess is required. 

V.A.3.b. A Technique to Construct an Initial Trajectory Guess 

An air-breathing SSTO HSV is likely to have three distinct flight phases when subjected to 
a dynamic pressure upper bound. The first phase seeks efficient performance in the subsonic 
to low supersonic flight regime while remaining below the dynamic pressure bound. When the 
bound is reached, the vehicle is likely to fly along it while accelerating to high hypersonic velocities. 
Finally, there is a climb to orbit, utilizing the rocket as necessary, while experiencing declininq 
dynamic pressure. An appropriate procedure to construct an initial trajectory guess is to develop 
feasible trajectory segments in each of the three identified flight phases An engineer’s 
experience-based intuition can be better applied to determination of likely vehicle flight 
characteristics in the simplified envelope of each of the separate flight segments. Solution of the 
segments in sequence from the ground up assures matching conditions between each seqment 
and provides a suitable initial guess for the overall trajectory analysis problem. 

A separate, simplified, optimization problem has been posed for each flight segment In 
each case, the propellant mass consumption has still been minimized. However, it remained 
necessary to construct an initial control history guess on each segment which flies a trajectory 
that nearly satisfies the locally applicable boundary conditions (though without reqard to 
performance). To simplify this development, the following assumptions have been made: 

- The active control used for pitch attitude is a rather than 6 e and 6 K . Elevon aerodynamic 
effects are retained by requiring maintenance of moment trim. Use of a (combined with 

<f> a and <P g ) is more intuitive for raw trajectory construction than 6 , and 6 K histories. The 
resulting control vector is: 
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- Rotation dynamics are not tracked (reducing the state dimensionality). The resulting state 
vector is: 


f 



x = 


I' 



(Eq. V.12) 


- Only air-breathing propulsion 


is assumed available on the first two segments, eliminating 

•» t- X*- >i> n »-> i ilr»l» m nhnoinn onn 


V.B. Demonstration Case Preliminaries 

V.B.I. Vehicle/Environment Models 

A 4200+ line file of data was used [11] to define the aerodynamic and propulsive char- 
acteristics, as well as mass properties of the example winged-cone hypersonic vehicle (HSV). 
Aerodynamic effects of rotational control aerosurfaces are included. The data is in table look-up 
form, involving up to trivariant dependencies for specified parameters, including some which 
evidence highly nonlinear behavior with respect to one or more variables. 

V.B.I. a. Available Database 

Variable Definitions and Dependencies: The variables in tables V.1 ,2,3 are included in the example 
vehicle model: 


Table V.1 

Example Model Aerodynamic Force Coefficients 


Variable 

Functional 

Dependencies 

Description 

C»u 

At , a. 

Drag increment coefficient for the basic vehicle 

(- D . <i « 

At . a , 6 

Drag increment coefficient for the right elevon 

(- D . tie 

At . a , 6 

Drag increment coefficient for the left elevon 

£ D .dr 

At , a , 6 

Drag increment coefficient for the rudder 

(-D.de 

At , a , 6 

Drag increment coefficient for the canard 

C,. a 

Af . a 

Lift increment coefficient for the basic vehicle 

(- L.d a 

At , a , 6 

Lift increment coefficient for the right elevon 

(-L.de 

At, a , 6 

Lift increment coefficient for the left elevon 

(-l.de 

M . a , 6 

Lift increment coefficient for the canard 

r . 

At . a . 6 

Side force increment for the riaht elevon 


































Table V.2 

Example Model Aerodynamic Moment Coefficients 


Variable 

Functional 

Dependencies 

Description 

r 

l . it a 

M . or . f> 

Rolling moment coefficient for the right elevon 

l . <t f> 

A / , a , t) 

Rolling moment coefficient for the left elevon 

l, dr 

M , a ,(j 

Rolling moment coefficient for the rudder 

c 

n ( n 

1/ , a 

Pitch moment increment for the basic vehicle 

r 

m . d a 

A/ , a . 6 

Pitch moment increment for the right elevon 

r 

/« . d e 

A7 , a , 6 

Pitch moment increment for the left elevon 

r 

v in , d r 

A / . a , fi 

Pitch moment increment for the rudder 

in . dr 

A/ , a.f) 

Pitch moment increment for the canard 

c,„ q 

M , a 

Pitching moment with pitch rate derivative 

r 

n . d a 

a; , a , 6 

Yawing moment increment coefficient for the right 
elevon 

r 

n.de 

A / , or , 6 

Yawing moment increment coefficient for the left 
elevon 

C n . ltr 

A / , a , 6 

Yawing moment increment coefficient for the rudder 

Cn <> 

A/ , a 

Yawing moment with roll rate dynamic derivative 


Table V.3 

Example Model Propulsion and Mass Property Coefficients 


Variable 

Functional 

Dependencies 

Description 

C T 

M . <|> , Q 

Thrust coefficient 

'.p 

Af , <t> , g 

Specific impulse 

Y 

eg 

m 

Longitudinal distance from moment reference center 
to vehicle center of gravity 

I yy 

m 

Pitch moment of inertia 
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Observations About Supplied Data: A number of parameter identities apply to the example vehicle 
model as a result of some of the symmetry used in its formulation: 


C D .da 
r 

rn . da 


C D . de 

r 

^ m . d e 


C / ,d(i C l, de 

c Y .da = ~Cy.de 


C L, da C 1 1 cfc 
C n , d a “ 


C /), dr(^) ^£>,dr( 

Some of the variables are nearly perfectly linear with respect to variations in angle of attack 
(a), eliminating the need to fit the data with respect to this parameter. These variables include: 

C y t d f . C l, dr ' C m . d r * ^ n , <i r * C n p 


The canard is assumed to be retracted above Mach 1 , so data are only provided for its 
aerodynamic contributions under subsonic flight conditions. 

Environment Model - Atmosphere: A standard 1976 atmosphere model that only treats attitude 
dependencies for physical properties has been used. Data are tabulated at one kilometer altitude 
increments. An exponential fit between altitude break points has been used. The atmosphere 
is assumed to rotate uniformly with the Earth’s surface directly below a specified location. 

Environment Model - Gravity: The gravity is modelled as an inverse square field based on the 
assumption of a homogeneous spherical Earth. Oblateness effects were ignored. 


V.B.I.b. Simplifications 

Reduced Rotational Dynamics: The trajectory analysis under this task was restricted to include 
the pitch plane rotational dynamics in order to reduce the computational burden required to 
accomplish the G&C analysis of the NASA/LaRC HSV design concept. This amounts to con- 
sideration of constant heading angle cases only, ignoring side forces. With this assumption pnly 
the following variables are retained, where left and right elevon contributions are merged into 
combined coefficients (by adding them together) since they act in unison for pitch control: 

C [>a > (C D. da + C D . de) • C D . dc ' C La' (^- L. da + C L.d s') ' C L.dc < 

C m a' (.C m.da + G mda'ifC mdc .C I sp , I yy , X cg 

Consistent with the exclusion of vehicle roll and yaw dynamics, the pitch dynamics are 
treated in a strictly planer sense. The vehicle wing axis is assumed to always be level with the 
local horizontal, and no sideslip is treated. Inertial coupling and frame rotation effects can 
eventually be added when 3-dimensional rotation dynamics are included (along with the 
aerosurfaces to provide the other rotational degrees of freedom). 

Exclude Canard Use: The canard in the example vehicle model is assumed to be usable only at 
low velocities (data ends at Mach = 0.9). It is excluded in the analysis methodology demonstration 
cases. All pitch control is accomplished with elevon usage. 

Curve Fits for the Vehicle Mass Properties Model: Since the data for the moments of inertia values 
and center of mass location included in the example winged-cone vehicle model was quite smooth, 
it was decided to develop a polynomial fit of their functional dependence. 

For the relationship of x cc to mass, a quadratic fit gave good accuracy. The resulting 

function is: 

x cc - -24.227373830 + .00008264665m - .0000000001 1 m 2 (^q. V.14) 

where mass is expressed in pounds and x cc in feet. In the analysis code, the masses are in 
slugs. Figure V-2 shows the fit curve and the data points (noted by circles). 
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A cubic fit was chosen for the / yy data, with the following result: 

/ yy = 2.8896578089m + ,00028207479m 2 - .00000000061m 3 


(Eq. V.15 ) 


The fit function for / yy is shown in figure V-3. 

Propulsion Phasing: The example vehicle data provides an air-breathing propulsion model to 
orbital velocities. However, to achieve the desired 110 NM circular orbit, some propulsive capability 
is required above the sensible atmosphere. A constant specific impulse rocket with throttle and 
a maximum thrust magnitude is assumed. Separate, discrete air-breathing engine off-time and 
rocket on-time values are defined as optimization switch points. The rocket is assumed to stay 
on from computed on-time to final orbit insertion. If a split rocket thrust history is necessary 
(similar to a Hohmann transfer), then intermediate zero throttle settings will result from the 
optimization solution. 

Rotation Effector Phasing: Use of aerosurface effectors for pitch control will be terminated at 
rocket ignition with control moments subsequently provided by rocket gimbal rotation of the thrust 
vector. Changes in rocket thrust along the vehicle longitudinal body axis will be ignored and 
normal force components due to gimbal rotation will not be included. 

Unconstrained Orbit Phasing: The target circular, polar orbit will be assumed to have an 
unconstrained phase angle. Any circular orbit of proper altitude will be deemed an acceptable 
result. 


V.B.2. Strategy for Final Near-Fuel-Optimal Trajectory Development 


V.B.2.a. Trajectory Definition (Boundary Conditions) 

Initial and final state conditions need to be set to complete the specification of the desired 
optimization problem. The initial condition is defined just after vehicle takeoff: 


*( 0 )- 



Ho 


Go 



(Eq. V.16) 


Only three of the terminal states are constrained as the problem is posed (including 
simplifications). The target orbit velocity is used to determine the forward integration cutoff 
condition Cl . The target orbit altitude and a zero target vertical velocity (for a circular orbit) fill 
out the equality constraint set, defining the two element constraint vector V . 


Q = ( V \ ) 2 + ( v 2 + 6 (AO e COS p.) 2 - (V targ ) 2 


(Eq. V.17) 


21 



(Eq. V.18) 


V = 


^ / ( l ) f targ 

^ l h ( x ) y 


V.B.2.b. Nine State Vehicle Dynamics Model 

A nine state vehicle dynamics model permits treatment of longitudinal rotational dynamics 
usina the aerosurfaces (elevons in the example vehicle) for rotational control. TTie state equations 
consistent with the LVLH coordinate system shown in figure V.4 and the modelling simplifications 
noted in section V.B.1 .b can be written as follows: 


/ 


- -V* 
■Vs 


\ 




v i ros.V' 


,v 4 

V, 


v 4 .v 6 - ( V S ) 2 tan x 


? F > 
a> „ v m v , - v , ( a> „ ) si n x 3 cos x 3 + — 

* 0 


x j- - - - - - — — + 2a> a (x 4 sin x 3 + x 6 cosx 3 ) + ~ 


X 


( ( v 4 ) +(X;i — -2u>,x s COS X - X , ( U) . cos x 3 ) 2 *~r 


\ 


M, 


i ( Q c t + r k ‘f* r 
f/o W *P« 1 sp » . 


(Eq. V.19) 


where: 



(Eq. V.20) 
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Using the zero bank angle assumption, the following force equations in the vehicle LVLH frame 
result assuming an arc tangent function definition over the range from -90 deg to +90 deg: 

f F i\ f )cosip(/ r 7 -cosY+ F N s\n y) 


F = 



si(j/t(v | )sin ip ( F r cos y + F w sin v) 
V -F r r,\n y+ F N cosy + F g 


(Eq. V.21) 


where sign ( v , ) has a value of 1 for v , > 0, and a value of -1 for v , < 0 . Also, note that only two 
force components remain in the body frame (no side force): 


F , = /’cos a - D 
F N = -(V sin a + /.) 

The following equations are used to relate flight angles to states: 

a = 0 - y 

.( 


Y = tan 


- A 


V V ( v 4 ) + ( x s)‘ 
v 


, I / ' v 5 

ij; = tan I — 


(Eq. V.22) 
(Eq. V.23) 

(Eq. V.24) 
(Eq. V.25) 

(Eq. V.26) 


The lift, drag, pitch moments, and engine thrust can be written in the following form involving 
coefficient functions specified in the example model data tables: 


where: 


l- = C L QS rt/ 

(Eq. V.27) 

D-C 0 QS r „ 

(Eq. V.28) 

C m Q c S rBf - .v cc ( /J>sin a + Lc os a) 
+ ('*-.v cc )<t> f T t sin6, 

(Eq. V.29) 

T-QC t + *,T, 

(Eq. V.30) 

c L = c La + c L ' 

(Eq. V.31) 

C D = D a + C l \ 

(Eq. V.32) 


(Eq. V.33) 
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(Eq. V.34) 


Uf = ^(^,) X + (^2) 2 + (^ 3 ) 2 

When using the example model database, v, and a are required, where v, is used along 
with r in an atmosphere model to compute M , . 

It is intended that eventually a model of moments from air-breathing engine thrust as a 
function of Mach number and angle of attack will be developed, the use of which will require 
additional terms in the vehicle pitch moment equation above. 

The gravity force acting on the vehicle can be written in the following form: 

= t l e - v 9 (Eq. v.35) 


The following equations are used to transform the throttle controls on the propulsion to 
unbounded control vector components. They are based on the assumption that < f > a ^0 and 

0 < <£, < 1 . 

4>a = ( lL a) 2 (Eq.V.36) 

<0,- - 5 {Eq. V.37) 

1 + 


V.B.2.C. Optimization Algorithm Formulation for the Example Application 

To perform closed-loop G&C studies on an air-breathing HSV during powered ascent to 
orbit it is necessary first to construct a representative, realizable trajectory under nominal flight 
conditions. The most appropriate way to construct the desired trajectory is to perform an 
open-loop analysis with the help of a trajectory optimization algorithm. Since the flight path must 
accommodate physical constraints such as limits on dynamic pressure, the methodology must 
directly treat the constraints to obtain a realistic trajectory. Also, since discrete phasing of pro- 
pulsive components is likely, the trajectory generation methodology must treat discontinuities in 
the system dynamics. A steepest descent, gradient-based optimization algorithm that treats all 
these concerns, as well as other features that would allow vehicle configuration characterization 
within the optimization, has already been developed and demonstrated [2]. To use this technique 
in the example application, it is necessary to formulate the optimization algorithm in an applicable 
form, and to define the vehicle dynamics in a compatible state space representation. 


The Algorithms Used: The trajectory optimization algorithm is based on a gradient/steepest 
descent technique for solving two-point boundary value problems. Predicated on the techniques 
reviewed in greater detail in refs. [2,4,12], the algorithm can accommodate equality and inequality 
constraints while optimizing performance with respect to state. The method can also be gen- 
eralized to treat variations in vehicle design, also as discussed in refs. [2,4,12]. To address the 
optimization problem, the vehicle state equation is first defined for use in the algorithm forward 

integration: 


clx 

— - = /(.Y,U,f s ) 

dt 


(Eq. V.38) 
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Then the algorithm development starts by definition of an applicable performance index: 

J ( t ) — <t> ( -V ( 'i: ) ) + f L(x(t),u(t))dt (Eq.V.39) 

J o 

A forward integration cutoff condition is defined (e.g. target orbital velocity): 

.Q ( v ( C ) ) = 0 (Eq. V.40) 

A Hamiltonian is constructed after defining a costate vector X which is determined by 
backwards integration: 


II = L+ \ T f 


(Eq. V.41) 


where: 

(It 

and where the backwards integration boundary condition is evaluated at the forward integration 
cutoff time x : 


X\x) 


4 >* 




V n x f J 




X(X),U(X> 


(Eq. V.43) 


An equality constraint vector is constructed to accommodate the terminal state boundary 
conditions defined at the forward integration cutoff time x : 

M / (.v(t)) = 0 (Eq.V.44) 


An additional adjoint variable, also determined by backwards integration, is defined to factor 
in the equality constraint functions (boundary conditions) at the terminal time: 

dA T (Eq. V.45) 

= “ f x A 

dt 1 x 

where the backwards integration boundary condition is evaluated at the forward integration cutoff 
time x : 


A r (x) 



n x / J 


X ( X ) . U ( X ) 


(Eq. V.46) 


It is helpful to introduce the following variables for subsequent notational simplicity where 
U and V are diagonal matrices of coefficients used to weight the relative size of perturbations to 
the controls and propulsive phasing switch times, respectively: 
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(Eq. V.47) 


M 


f ^ ( * s 1 ) ( / s l f S ! ) \ 

f's 2 ) 

0 = 

V : J 

a 7 (( s ,)(/;,- A s ',)i-A r (t, 2 )(/; 2 -/; 2 )i-] 
l„ y - r^ T f u u il/\dt 

J 0 

f X A T f u UHldt 
J o 


(Eq. V.48) 


(Eq. V.49) 


(Eq. V.50) 


ijj- \ x y\ u uH[dt 

J 0 

/i = (l tv + M/M 7 )' 1 

c - (i VJ + 

The resulting algorithm variational equations for the control vector u and the 
phasing switch times t s (where dynamic discontinuities occur) are: 

bu = -U{f T u /\BC^ + C^U T a - f T u J\BC)) 
bt s = -V{M r BC^V + C S (g-M T BC)) 


(Eq. V.51) 
(Eq. V.52) 
(Eq. V.53) 
propulsive 

(Eq. V.54) 
(Eq. V.55) 


where: 

dJ s + (\ T vjB + g T VM T B)C W V 
^ J = ~ \ jj - \ y } BC + g T Vg - g T V\A T BC 


(Eq. V.56) 


With specification of ci J s and C the variational equations can be used to perturb the 
desired control history on successive iterations through use of Eqs. V.54,55 to eventually achieve 
vehicle performance optimality while satisfying all applicable constraints. The scalar dJ s reflects 
the specified cost improvement sought on each optimization algorithm iteration. The diagonal 
matrix C ^ weights the relative size of equality constraint improvements on each optimization 

algorithm iteration. Proper selection of these weightings as well as U and V is critical to the 
numerical stability and convergence rate of the optimization process. A block diagram showing 
the signal flow and processing order of the optimization algorithm is shown in figure V-5. 
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Treating Constraints: There are several types of constraints applicable to HSV trajectory opti- 
mization. As is true for any two-point boundary value problem, it is necessary to apply equality 
constraints at the trajectory final target condition. These are treated by adjoining the variable A , 
thereby incorporating the influence of the boundary conditions. 

There are constraints on some of the controls. An example would be a lower bound (fuel 
flow off) on throttle setting. Often such constraints can be treated by nonlinear mapping into a 
control space without bounds. Care must be taken when using nonlinear mapping to avoid initial 
solution guesses that specify mapped control values in insensitive regions of the actual control 


Inequality constraints need to be applied as a result of physical design constraints such as 
dynamic pressure and thermal limits. For this methodology, penalty functions are introduced 
through the integral performance index L . This approach requires periodic adjustment of penalty 
function gains as the optimal solution is approached. 


Required Special Features: As is generally true when two-point boundary probjem optimization 
a gorithms are used, some application specific considerations are required to assure qood 
algorithm performance for the HSV problem. These include the careful choice of the integration 
cutoff condition as well as special approaches to balancing the perturbations to boundary value 
constraint violations against perturbations to controls needed to improve performance. Failure 
to carefully treat these issues can prevent acceptable algorithm numerical behavior. 


Unlike a rocket ascent to orbit, HSVs are more likely to follow a trajectory at low Mach 
numbers akin to supersonic fighter aircraft minimum time to climb profiles. Often this results in 
a zoom maneuver" in the transonic range, where the vehicle uses already acquired potential 
® ner 9y 1® coml pined with thrust, to quickly transit the maximum drag environment near 

Mach 1 [14], The resulting ascent can involve a negative vertical velocity component for a short 
interval me nonmonotonic behavior of this velocity component makes it an unacceptable 
parameter for algorithm integration cutoff. Either horizontal velocity or orbital energy state provides 
much greater likelihood of achieving acceptable state conditions at integration cutoff. 

• Successful convergence of the optimization algorithm to the desired minimum cost solution 
requires that the effort to drive the violation of the boundary conditions to zero be balanced aqainst 
the effort to find an optimal control history. Each parameter requiring perturbation must be 
separately weighted, with the emphasis that each receives chosen to stress those variables that 
currently promise the most improvement without numerically destabilizing the alqorithm. (Bad 
choices could induce highly nonlinear steps.) As the interim trajectory solution approaches 
optimahty, the relative sensitivity of the perturbed parameters change. Techniques must be 
applied to determine current sensitivity and modify the weightings as the algorithm progresses 
In some instances, this may require acceptance of perturbations that do not uniformly improve 
all boundary value violations. A composite measure of proposed algorithm step acceptability 
can be constructed that takes precedence over individual variation effects, as discussed in ref. 


Cost Function Specification: Of primary importance in applying an optimization algorithm to a 
trajectory analysis problem is specification of the terms of a cost index. For this task the terminal 
cost term was used to minimize vehicle propellant mass usage to go from the ground to orbit 
and the integral cost term was used to constrain the dynamic pressure. The resulting functional 


<j> = m 0 -m(T) (Eq. V.57) 

L = C q (Q-Q d ) 2 u 0 [Q-Q d ] (Eq.V.58) 

where « 0 fQ ~ Qd] is a unit step function that has zero value for Q <Q D and has a value of 1 
otherwise. 
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V.B.2.d. Relevant Partial Derivatives 

Boundary Condition Derivatives: The following two equations are the example specific derivatives 
of the boundary condition equality constraint functions needed in Eqs V.43,46. 

2 uo e cos 11(^2 + rco e cosp) \ 


f 


0 


n T x = 


\ 


- 2roo e sin \i(u 2 + roo e cos \i) 

2t/, 

?.(v 2 + roo e cosp) 

0 
0 
0 
0 


(Eq. V.59) 


_f\ 0000000 0 ' 

0 0 0 0 I 0 0 0. 


(Eq. V.60) 


Cost Function Derivatives: The following two equations are the example specific derivatives of 
the performance index function J needed in Eqs. V.42,43: 

<j) x = (0 0 0 0 0 0 0 0 -1) 


dQ 

L v = 2.C q(Q ~ Q q[Q ~ Q d] ^ x 


(Eq. V.61) 
(Eq. V.62) 


where: 


d -Q= ( 0 0 pu, p o 2 P o 3 0 0 o] 

dx \dr 2 ' 


(Eq. V.63) 


and where is assumed to be a constant. The values for p and £must be derived from a 
standard atmosphere model. 

The Hamiltonian derivative with respect to controls needed in Eqs V.50,51 ,54 can be 
expressed by the relationship H u = L u + ^ r / u which requires a partial of the integral cost function 
with respect to the control vector. For the function of L in the example here, L u has zero value. 

State Equation Derivatives: The partials of the state equation f in Eq. V.19 with respect to states 
x is required for use in Eqs. V.42,45. Because of the number of equations required to complete 
the definition of / x for the example model, details are provided in appendix A.1 . 

The partials of f with respect to controls u is required in Eqs. V.49,50,54. Because of the 
number of equations required to complete the definition of / u for the example model, details are 

provided in appendix A.2. 
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V.B.3. Preliminary Near-Fuel-Optimal Trajectory Development 

A general methodology to perform integrated trajectory and control analysis on SSTO HSVs 
has been devised, and the mathematical basis for demonstrating the complete minimum-fuel- 
to-orbit trajectory analysis of the example vehicle has been presented above. However, to use 
the methodology, it is first necessary to construct a complete, mathematically realizable trajectory, 
from horizontal takeoff to orbit, which initializes the optimization algorithm. While manual iteration 
of the vehicle dynamics model could be performed to provide an overall trajectory guess, the job 
can be made much less labor intensive by breaking the problem into a series of simpler initial 
trajectory development problems (each still requiring a manually developed starting trajectory 
profile). Sequential development of the trajectory segment guesses from the ground up supports 
local optimization of each segment, resulting in a composite trajectory that is likely to be better 
behaved when used to initiate the entire ground to orbit trajectory optimization analysis. 

The formulation of three representative trajectory segments are presented here to support 
the initial complete trajectory definition demonstration. To further reduce the computational 
complexity of the initial example trajectory segment development, rotational states of the vehicle 
are not tracked, and the moments are only included to treat vehicle lift and drag effects which 
result from maintaining pitch trim throughout flight. Engine throttle and vehicle angle of attack 
are treated as the controls (not aerosurface and rocket gimbal deflection angles). A simpler 
seven state vehicle dynamics model results. 

V.B.3.a. A Representative Trajectory Partition 

Three ascent trajectory legs were defined, each with special characteristics, which simplify 
development of an initial trajectory guess. The first leg flies from a horizontal takeoff up until the 
dynamic pressure constraint is reached at a target velocity in the Mach 2 range. The second leg 
flies as close to constant dynamic pressure as possible. The third leg ascends to the orbital 
condition from the end of the constant dynamic pressure segment. The following subsections 
further discuss each leg. 

Leg 1 - Takeoff to the dynamic pressure limit: This segment is solved first using the desired 
vehicle horizontal takeoff conditions as initial conditions. A target terminal time velocity is defined, 
selected on the basis of the desired dynamic pressure limit (near Mach 2 for a dynamic pressure 
bound of 1 000 psf). Terminal time constraints on the trajectory segment conditions are selected 
so the dynamic pressure bound is approached with near zero slope to provide proper initial 
conditions for the next leg which seeks constant dynamic pressure. Inequality constraints are 
applied only if required to keep the peak dynamic pressure close to the terminal constraint value. 
Propulsion is limited to the air-breathing system. This results in the elimination of the switch points 
in the trajectory leg analysis. When applying the optimization algorithm g and M have zero value, 
and the variations 6 t s disappear. Also the following simplified equations result: 


*-( Iv.v ,)' 1 

(Eq. V.64) 

6ix--f/(/lAflC V( V I +C Ji (Hl-/lABI ViJi )) 

(Eq. V.65) 

C* i 

(Eq. V.66) 

lj j - Ij j B\y j 

JjJj M J j Y 1 J 1 



Leg 2 - Constant dynamic pressure: This segment is solved second, using the terminal conditions 
from leg 1 as initial conditions. A target terminal time velocity is defined, selected on the basis of 
where constant dynamic pressure air-breathing flight, with net acceleration, is still barely 
sustainable (assumed to be 20,000 ft/s air speed). The same terminal time constraints are applied 
as those used on the first trajectory leg. A distributed cost term is applied to the trajectory segment 
to nearly maintain the desired constant dynamic pressure for the entire leg. Propulsion is limited 
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to the air-breathing system. As for leg 1 , this results in the elimination of the switch points in the 
trajectory leg analysis. When applying the optimization algorithm A, g, M and the variations 6 t s 
disappear. Also the following simplified equations result: 


«-( 1» 2 » 2 )" 

(Eq. V.67) 

5u - -£/( f l A BC t2 M< + C j2 (Hl - fl A 

(Eq. V.68) 

</j s + i ; j BC t V 2 

(" 

2 lj j ~ I't j 

J 2 J ?, T 2 J 2 r 2 J 2 

(Eq. V.69) 


Leg 3 - Final Ascent to orbit: This segment is solved third, using the terminal conditions from leg 
2 as initial conditions. The target circular orbit velocity is used as the terminal time cutoff condition 
with the orbit altitude and radial velocity also constrained at terminal time. The inequality constraint 
on dynamic pressure is still enforced as required. Propulsion can include the air-breathing system 
and the rocket. The mathematical form of the problem is the same as discussed above for the 
entire trajectory solution, but using W 3 as the terminal constraint function and J 3 as the cost. 


V.B.3.b. Simpler Seven State Vehicle Dynamics Model 

The simpler seven state vehicle dynamics can be written as follows in state equation form 
consistent with the frame of reference given in figure V.4: 
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where: 




X = 


^2 

V 


(Eq. V.71) 


v::y 


and where f R has zero value for ascent legs 1 and 2. For the orbit insertion leg 3: 




T k * k 


(Eq. V.72) 




Using the zero bank angle assumption, the following force equations in the vehicle LVLH frame 
result assuming an arc tangent definition over the range from -90 deg to +90 deg: 



( F 


F = 


= 




f sign(v , ) cos ip(F r cosy + f^sin y) 
sign (u,)simp( F r cos Y+^sinY) 
- F T sin y + F n cosy + F g 


(Eq. V.73) 


where sign{ v i ) has a value of 1 for v [ ^ 0 , and a value of -1 for v 
force components remain in the body frame (no side force): 

F T - T cos a - D 


< 0 . Also note that only two 
(Eq. V.74) 


F N - -(Tsin a + L) 


(Eq. V.75) 


Also, the following equations are used to relate flight angles to states: 


Y ** tan 




n /(* 4 ) 2 + Os )‘ 


(Eq. V.76) 


ip - tan 



(Eq. V.77) 


The lift, drag, pitch moments, and engine thrust can be written in the foljowing form involving 
coefficient functions specified in the example vehicle model data tables, utilizing the assumption 
that the elevons are used only to maintain a vehicle pitch trim condition. With moment balance, 
pitch rates must be assumed zero, eliminating the need for a moment coefficient due to pitch 
rates C m . This further simplifies the dynamics model considerably, but introduces some 

mathematical inconsistency since a time variant angle of attack actually implies vehicle rotation 
rates: 
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(Eq. V.78) 
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ef 


where: 


D-C d QS„, 

(Eq. V.79) 

= C fll QrS n>f - .v cc (Dsin a + L cosa) 

(Eq. V.80) 

T-QC r + * K T K 

(Eq. V.81) 

Ci- C L + C l 

(l 0 

(Eq. V.82) 

C d = c D + C D 

a q 

(Eq. V.83) 

C m = ^„ + C me 

(Eq. V.84) 

V( i 'l) 2 + («2) 2 + (".l) 2 

(Eq. V.85) 


When using the example vehicle model database, M f , a, and 6 c are required, where v f is 
used along with r in an atmosphere model to compute M f . An implicit equation to determine 
6 e is derived from the moment balance relationship. Only the elevon drag, lift, and moment terms 
are dependent on 6 e . Other terms, and the other degrees of freedom for the elevon dependent 
terms, are determined by the currently stored history of a , m , and M f . Using the following 
equation (and the associated relations defining the equation coefficients), a value of 6 e can be 
found iteratively. 
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The gravity force acting on the vehicle can be written in the following form: 

= \i a x 7 
,fJ U ,) 2 


(Eq. V.86) 

(Eq. V.87) 
(Eq. V.88) 
(Eq. V.89) 

(Eq. V.90) 


Finally, the following equations are used to transform the throttle controls on the propulsion 
to unbounded control vector components. They are based on the assumption that <f> a >0 and 

0< <t> R < 1 . Also, 4> R only applies to leg 3 of the ascent trajectory. 
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(Eq. V.91) 


*„ = (' o 2 

( N = T ~ — 3 (Eq.V.92) 

1 + O/r) 


V.B.3.C. The Boundary Conditions 

Boundary conditions are set for each of the trajectory legs based on the criteria discussed 
above. Velocities are used to identify the end of each trajectory leg . Relative air speeds are 
used for legs 1 and 2, and Earth centered inertial velocity is used for leg 3 (assuming no radial 
velocity term). Other terminal time conditions (Y , ) apply to trajectories 1 and 3. For leg 1 , the 

dynamic pressure limit to be applied on leg 2 must be matched, and have zero slope. For leg 
3, the additional constraint conditions assure a circular orbit at the desired altitude. If a Hohmann 
type transfer is required on leg 3 (two separate burns with the later being a circularization impulse) , 
then the optimization solution will accomplish that result by forcing the rocket throttle setting to 
zero for a coast interval. 
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(Eq. V.93) 
(Eq. V.94) 
(Eq. V.95) 


(Eq. V.96) 


(Eq. V.97) 


(Eq. V.98) 


V.B.3.d. Cost Function Specification 

Of primary importance in applying an optimization algorithm to a trajectory analysis problem 
is specification of the terms of cost indices. The basic function forms are: 


J|(T,) = 4>,(v(x,))+ f L l (x(t),u(t))dt 

J 0 

(Eq. V.99) 
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(Eq. V.100) 
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J .( i ,)' «,( '(! ,)) - f J L 3 (x(t),u(l))dt 

J t 2 

where: 
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(Eq. V.101) 


(Eq. V.102) 


Note that a 2 is only used in leg 3 of the ascent trajectory. 

For the current example, the cost terms used are as follows: 
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(Eq. V.103) 
(Eq. V.104) 
(Eq. V.105) 

(Eq. V.106) 

(Eq. V.107) 
(Eq. V.108) 


V.B.3.e. Relevant Partial Derivatives 

Boundary Condition Derivatives: The following equations are the required derivatives of the 
boundary condition functions in Eqs. V.93-98 used in Eqs. V.43,46: 
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(Eq. V.109) 
(Eq. V.110) 


(Eq. V.111) 
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(Eq. V.112) 


where: 

dQ 

d v 




/ a££C|) \ 

ix 

iQX<) 




f i o 
^0 0 


/ <*<?( ll) \ 

<K?(T,) 

V / 

0 0 0 0 

ooo 1 


o 

o 


dpO/)^ 
dr 2 


0 0 p v , 


P V 2 P^3 



dQ 

dt 


f i / ) (tp , r r rN 
^ — + P(v 4 / 4 + - v 3/5 + x 6.fa) 


(Eq. V.113) 


(Eq. V.114) 


(Eq. V.115) 
(Eq. V.116) 


/ v 2 f d 2 p d p \ 


ff(I)Y 

^ dx 


\ 


2 dr 


0 

0 

dp 

-X 4 X,- + p/ 4 
dp 

-.V- r> .V 6 - + p/ 5 

_3 x j£ 5£ + p/4 

2 dr 6 
0 


f if* if s ^/ 6 V 

P X 4 + -V,; +X A 

V v 3x 6 dx J 


/ 


(Eq. V.117) 


where the last term of the previous equation was written in vector form for notational simplicity. 
Also, note that since p is assumed to be only a function of r , y r = j- r . 
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Cost Function Derivatives: The following equations are the derivatives of the cost functions J , 
defined in Eqs. V.103-108 that are required in Eqs. V.42,43: 
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(Eq. V.118) 
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(Eq. V.119) 
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(Eq. V.120) 

For the first segment, 1 . , v 

= 0 if Q < Q d , and when Q > Q D : 
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(Eq. V.121) 
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For the second segment: 
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(Eq. V.122) 


For the third segment, as for the first, L 3 v - 0 if Q < Q D , and when Q> Q D : 

30 (Eq. V.123) 

l. ;li = ZA (? 0 ) — 

The Hamiltonian derivative with respect to controls can be expressed by the relationship 
1 1 „ = L u + X r f u which requires a partial of the integral cost function with respect to the control 
vector. For the functions of I selected in the example here, L u has zero value. 

State Equation Derivatives: The partials of the seven state vehicle dynamics equation f given in 
Eq. V.70 with respect to states is required for use in Eqs. V.42,45. Because of the complexity of 
its specification, details are given in Appendix B.l. 

The partials of the seven state vehicle dynamics equation f given in Eq. V.70 with respect 
to controls is required for use in Eqs. V.49,50,54. Because of the complexity of its specification, 
details are given in Appendix B.2. 


V.B.3.f. Calculation of Elevon Trim Deflection 


Calculation Method for Elevon Trim Deflection: A general numerical method has been devised 
to calculate an aerosurface deflection that trims a vehicle’s moments resulting from aerodynamic 
loading. In the example application, the method has been adapted to calculate elevon trim 
deflections. The vehicle longitudinal moment balance equation is given by: 


QS r „)C{C m n + C m e ) = Q S rff x cg ((C D a 


C 0 e )sin(a) + (C La + C L Jcos(a)) 


(Eq. V.124) 


where the left side is the intrinsic aerodynamic moment, and the right side is the added moment 
caused by center of gravity offset. If x CG is zero, the trim deflection is found by obtaining the 

value of 6 a such that: 
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(Eq. V.125) 
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When x cu is nonzero, as is generally the case, the flap deflection enters the balance equation 
through C L a and C D a in addition to C m a . To find the value of 6 e that satisfies the above equation 
we rewrite the equation in the form: 

K,C B '-K 2 C L '-K 3 -0 (Eq.V.126) 

where the k , coefficients are independent of 6 e . We then apply a Newtonian iteration that finds 
a zero of the above expression for the solution value of 6 e . The iteration is of the form: 

(Eq. V.127) 
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(Eq. V.128) 
(Eq. V.129) 


To speed convergence, the previous solution for 6 C is used as the starting point for 
successive time steps, and the resulting convergence typically requires just a few iterations. To 
stay within the 20 degree limit on angle of attack, the starting perturbation values of 6, that are 

used to establish the derivative d f (6 e ) are set one degree towards zero from the previous trim 
position. For example, if the previous trim was -16 degrees, the perturbation point is chosen as 

~ 1 o ucyrsGS. 

Calculation Method for Elevon Deflection Partial Derivatives: Close examination of the elevon 
deflection partial equation reveals that on the right hand side, the partials — and — are 

d x d x 

themselves functions of ~ . To solve for 77 we write: 
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where i is an index over the states, and the coefficients C are given by: 
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(Eq. V.131) 
(Eq. V.132) 

(Eq. V.133) 
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Next we expand ~~~ and — in a similar manner, namely. 
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(Eq. V.134) 


(Eq. V.135) 
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(Eq. V.137) 
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(Eq. V.139) 
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Substituting the expansions for ^ and into the expansion for ix gives the equation. 
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which has the desired solution: 
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V.B.3.g. Initial Conditions 

The initial conditions assumed for the example vehicle in the analysis demonstration are 
set immediately after horizontal takeoff: 



20,909,723// 
280 dpg 
28.5225deg 
444.889 ft/s 
-37 .'5698 ft / s 
-6.46676 ft/s 
299,000/6 


(Eq. V.142) 
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VC. Analysis Demons tration Results 
V.C.I. Cases Run 

An analysis methodology demonstration case was run using the example vehicle model. 
The case was run from the moment of horizontal takeoff to 20,000 ft/s relative air speed. The 
traiectorv ootimization objective was to achieve a near minimum propellant consumption condition 
while applying a dynamic pressure constraint of 1 000 psf. The Presented results are an aggregate 
of computation on two flight legs: Takeoff to Mach 2 (leg 1); Mach 2 to 20,000 ft/sec (Ieg 2). The 
second leg was constructed from two shorter segments which were merged at 10,000 ft/s. 

Oriainallv the leg 1 results were run first, using a preliminary version of the example vehicle 
model database, and the terminal conditions for that trajectory segment were used as initial 
conditions for leg 2 analysis. After starting the leg 2 computations, the final version of the example 
vehicle model was prepared which amended the sign of the center of mass to moment reference 
center distance. Based on the model update, the leg 2 analysis was restarted, and run to 
completion. Subsequently an effort was also made, using the amended data set, to regenerate 
lea 1 results To get a reasonable new interface boundary condition between the two trajectory 
leas matches were mathematically forced on the dynamic pressure, its time derivative, and the 
vehicle mass at the interface point. (The dynamic pressure and derivative matching conditions 
had the effect of forcinq an altitude match given the combination of the atmospheric model used 
and the fixed trajectory segment cutoff velocity.) The effects of these matching conditions on 
trajectory results are discussed in section V.C.2.d. 

V.C.2. Observations from Demonstration Case Results 

Manv significant vehicle dynamics and control characteristics are apparent from inspection 
of the plotted data. These results are indicative of the benefits of the generic integrated 
trajectory/control analysis methodology. There are features that provide evidence of flight 
characteristics that should be treated explicitly rather than by being buried in data tables (par- 
ticularly discovered propulsive discontinuities). Some results provide information essential to 
closed-loop G&C system design requirements development. Also, vehicle performance 
sensitivities to modelling uncertainty have become apparent. 

V.C.2.a. Observations from Specific Variable Plots 
Angle of Attack 

Fiaures V-6 7 plot the computed angle of attack vs Mach number and time, respectively. 
The history shows the high angle of attack at takeoff necessary to get sufficienl .lift. It ® 
local minimum in the transonic region allowing the vehicle to use more of available engine thrust 
to accelerate through the high drag phase, rather than fighting lift related drag and gravity. 
Generally thereafter a gradual increase occurs to Mach 6. Above Mach 8 the angle of attack 
declines smoothly with Mach number as a result of orbital mechanics (the lift requirement declines 
as velocity approaches orbital speed). Some of the additional detail in the plot structure is 

discussed below. 

Equivalence Ratio 

Fiaures V-8,9 plot the engine throttle equivalence ratio vs Mach number and time, 
respectively. The history shows a local peak in the transonic region to get more thrust in the high 
draa fliaht phase. Another local peak occurs near Mach 6 where initial scramjet use is likely, and 
neither a ramjet nor a scramjet is capable of particularly good performance. Above Mach 8 the 
equivalence ratio undergoes very little variation. 
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Elevon Deflection 

Figures V- 10,1 1 plot elevon deflection vs Mach number and time, respectively. The history 
shows a maximum value near Mach 2 and reflects a change in vehicle trim as a result of combined 
changes in aerodynamic loading and center of mass shift due to fuel consumption. Above Mach 
8 elevon deflections eventually diminish as the angle of attack declines (figures V-6,7), offsetting 
the effects of center of mass shifts due to fuel consumption. 

Dynamic Pressure 

Figures V-1 2,1 3 plot the vehicle dynamic pressure and the dynamic pressure time derivative 
respectively vs Mach number. The profile shows a capture of the 1 000 psf constraint at about 
Mach 2. Most of the leg 1 flight phase is unaffected by the dynamic pressure constraint until very 
near Mach 2. By contrast, the leg 2 results reflect a near constant dynamic pressure flight condition 
objective. Note that the leg 2 results were still very slightly sub-optimal when computation was 
stopped, accounting for the ±2-4 psf oscillations in the dynamic pressure around the 1000 psf 
constraint. 

Flight Path Angle 

Figures V-1 4,1 5 plot the vehicle flight path angle and flight path angle time derivative 
respectively, vs Mach number. The flight path is very steep in the initial flight phases to gain 
altitude quickly, thereby assuring that the vehicle remains within the dynamic pressure constraint 
as it accelerates at higher Mach. As the vehicle accelerates toward orbital velocity, the angle 
tapers off to near zero. 

Vehicle Acceleration 

Figures V-1 6,1 7 plot the vehicle ground (also air) relative acceleration and inertial relative 
acceleration respectively vs Mach number. Both plots indicate that maximum acceleration occurs 
between Mach 2 and 6. Above Mach 8 the vehicle acceleration initially drops due to declining 
engine performance, but eventually starts rising again due to declining vehicle mass. While the 
inertial acceleration is much lower than the ground relative acceleration at low Mach (due to 
ground motion), both become almost identical at high Mach. 

Vehicle Mass Behavior 

Figures V-1 8, 19 plot the vehicle mass and mass flow respectively vs Mach number. The 
overall mass trend is a nearly linear decline with Mach number (as opposed to an exponential 
decay typical of a rocket). However, above Mach 8 the fuel mass flow magnitude goes through 
an oscillation that reflects the nonlinear propulsion system (specific impulse) model and 
associated curve fit. 

Altitude History 

Figure V-20 plots the vehicle altitude vs Mach number which shows the square root 
function-like shape characteristic of HSVs going to orbit. As the Mach number increases, less 
altitude rise is required to maintain the dynamic pressure constraint for a given velocity change. 

Mach History 

Figure V-21 plots the vehicle Mach number vs time. The acceleration drops near 400 
seconds at Mach 8 leading to a lower Mach increase slope for the remainder of the flight. 

V.C.2.b. Observed Acceleration and/or Mass Flow Discontinuities 

A number of discontinuities are observed in vehicle acceleration characteristics that 
probably reflect discrete transitions in air-breathing propulsive modes. The acceleration plots 
(figures V-1 6, 1 7) are consistent with the use of a multi-mode propulsion system, possibly a turbojet 
to near Mach 3.5, followed by a ramjet to near Mach 6, with a scramjet for the remainder of 
air-breathing flight. These characteristics were not explicitly identified before performing the 
numerical analysis. The strong nonlinearities that they induce caused a large increase in com- 
putational effort to achieve near-minimum-fuel trajectory convergence. If the propulsive phases 
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were treated explicitly in the analysis, there would have been less computational burden along 
with the added benefit of being able to assess the optimal timing of the propulsive phasing. Also, 
the separation of the specific impulse and thrust coefficient data in the example vehicle model 
makes possible inadvertent introduction of model inconsistencies by incorrect specification of 
individual data points. Some model inconsistencies do seem to have been identified by the 
numerical analysis of the trajectory. These unspecified discontinuities and modelling inconsis* 
tencies are evidence of a situation where there should be close coordination between developers 
of propulsion system models and G&C system analysts. Specific examples follow: 

- At Mach 3.5, there is a steep change in vehicle acceleration (see figures V-16,17) that 

correlates with a sudden change in fuel mass flow rate. Inspection of the example vehicle 
model propulsion data shows the cause of the odd behavior seen in the plots is a sudden 
change in modelled specific impulse behavior at Mach 3.5 (data is supplied at Mach 3, 
3.5, and 4). This may reflect a change from turbojet to ramjet operations. 

- A steady rise in the equivalence ratio begins at about Mach 5.45 (see figures V-8,9). The 

start of the equivalence ratio rise correlates with propulsion data that suddenly results in 
higher specific impulse values at equivalence ratios between 1 and 2.5 in the flight range 
between Mach 4 and 6. This may correlate with a transition from ramjet mode to a scramjet 
mode as Mach 6 is approached. 

- There are some significant data variations between Mach 7 and 8. Below Mach 7.5, a 

sharp local maximum occurs in the angle of attack (figures V-6,7), a rapid drop occurs 
in the equivalence ratio (figures V-8,9), and a local minimum occurs in the mass flow rate 
(figure V-1 9). At Mach 8.0, there is a nearly discrete change in vehicle acceleration (figures 
V-16,17) that correlates with a sudden change in fuel mass flow rate (figure V-1 9) and a 
nearly discrete change in the desired equivalence ratio (figures V-8,9). Inspection of the 
example vehicle propulsion model shows major changes between Mach 6 and 8 in the 
character of both the thrust coefficient and specific impulse tabulated data, somewhat 
akin to the behavior observed at Mach 3.5. This probably reflects some propulsion model 
fidelity problems and possible inconsistencies in the specific impulse and thrust coefficient 
data trends. While spline routines have been used to smooth out corners in the tabulated 
vehicle model data, the routines were not intended to effectively handle nearly discrete, 
multidimensional changes in functional characteristics. 

- The most significant remaining variations from the 1000 psf dynamic pressure constraint 

(figures V-1 2 13) are at Mach 3.5 and near Mach 6. These points correspond exactly 
with propulsion model data points and behavioral changes intrinsic in the tabulated 
example vehicle model data. Transients in the elevon trim angle (figures V-1 0,11), 
equivalence ratio (figures V-8,9), angle of attack (figures V-6,7), acceleration (fiQures 
V-1 6 17) and engine mass flow (figure V-1 9) occur at the same Mach numbers. This is 
more evidence of unspecified discrete changes in the propulsion model. The ground 
relative acceleration plot shows classic turbojet performance below Mach 3.5, with a 
sudden jump in acceleration (figures V-16,17) to a ramjet like profile. Above Mach 5.5, 
the acceleration profile looks more like a scramjet, peaking just beyond Mach 6. 


V.C.2.C. Analysis Observations Affecting G&C Design 

The power of the integrated trajectory/control analysis methodology to identify vehicle flight 
characteristics with critical impact on the desired closed-loop G&C system design is evident in 
some of the demonstration case data. Examples include large elevon trim deflections which 
indicate important limits in available effector control moments for broad flight regimes, and epi- 
sodic tendencies to seek rapid angle of attack variations which are evidence of a requirement for 
a hiqhly responsive controller. These two effects actually provide conflicting G&C development 
guidelines, creating narrow constraints on an acceptable design. Some specific observations 
based on the plotted data follow. 
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- The elevon deflection angle history (figures V-10,11) reflects changing vehicle trim char- 

acteristics due to Mach number dependent aerodynamic effects and vehicle center of 
mass changes that result from fuel consumption. The direction of vehicle trim moments 
reverses at about Mach 4. As the Mach number increases beyond Mach 4, the 
increasingly negative trim angles induce added drag that must be overcome by pro- 
pulsion settings. Much of the net vehicle drag losses at high Mach is due to the elevon 
trim deflection. Between Mach 1 0 and 1 5, the elevon trim angle exceeds 1 6 deg (peaking 
at 1 7 deg) leaving less than 4 deg authority for active control. 

- A steep perturbation in the angle of attack (figures V-6,7) occurs at about Mach 5.24, and 

reverses at about Mach 5.45. It seems to correlate with some large, highly nonlinear 
changes in the engine model specific impulse data that occur between Mach 4 and 6 for 
the equivalence ratio range (figures V-8,9) applicable to the profile. 

- The gamma dot plot (figure V-15) has steep variations between Mach 5 and 6 in direct 

correlation with the angle of attack variations (figures V-6,7) noted above. These variations 
have similar but smaller effect on vehicle acceleration (figures V-16,17) and mass flow 
(figure V- 19). 

- A sharp local minimum occurs in the angle of attack (figures 6,7) at Mach 6.35. It is also 

seen in the flight path angle time derivative data (figure V-1 5). The physical (model based) 
origins of this effect are not yet fully understood. 

V.C.2.d. The Mach 2 Transients 

As was noted above, an amendment to the vehicle model as its final version was prepared 
resulted in the assessment of the takeoff to Mach 2 leg and the Mach 2 to 20,000 ft/sec leg to be 
evaluated in reverse of the usual order for trajectory segment generation. As a result, and in spite 
of enforcement of several interface matching conditions, not all the state function and control 
histories were continuous across the interface. This behavior is evident in transients in some 
plots at Mach 2 that are identified below. The transients are confined to state and/or control 
functions not forced to satisfy matching conditions at the trajectory leg interface. If the very brief 
Mach 2 transients are ignored, the plots still give a good indication of the constrained near- 
fuel-optimal trajectory and control flight characteristics for the example vehicle. The specific Mach 
2 data problems seen in the figures are as follows: 

- There is a sharp angle of attack drop and recovery (figures V-6,7). 

- There is a sharp equivalence ratio drop and recovery (figures V-8,9). 

- There is a very brief drop and recovery in the elevon deflection angle (figures V-10,1 1). 

- The flight path angle shows a negative spike (figure V-1 4). With a full state match at the 

trajectory segment interface, this state function transient would not be physically possible. 

- The flight path angle time derivative (figure V-15) shows a spike corresponding the path 

angle transient. 

- There are sharp drops and recoveries of the vehicle acceleration (figures V-16,17). 

- There is a negative then positive change in the propellant mass flow rate (figure V-1 9). 

Results processed before the model amendment suggest that the positive direction 
transient is real. (It may be the result of a propulsive mode change.) However, the spike 
in the negative direction is probably due to the trajectory segment interface matching 
problems. 

V.C.2.e. Performance Sensitivity to Modelling Uncertainty 

Significant differences in the effective dynamics of specific HSV models can result from 
alternative methods of computing data table based performance. Analysis of the example vehicle 
model propulsion system data indicated that steep changes in response occur near specific 
Mach numbers. These changes can induce significant differences in engine performance 
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estimates in the vicinity of those data points when different interpolation schemes are used. 
Trajectories computed with spline fitted vehicle model representations were compared with tra- 
jectories computed with linearly interpolated models. One comparison showed nearly a 5,000 lb 
variation in propellant used between Mach 2 and Mach 3.5 as a result of the dissimilar interpolation 
strategies. 

The differences in vehicle performance estimates due to interpolation methodology have 
been shown to have the potential of drastically affecting the predicted vehicle payload to orbit. 
This indicates that the integrated trajectory/control analysis methodology can identify vehicle 
model fidelity inadequacies in sensitive flight regimes. It also is evidence that the methodology 
can determine when the specific vehicle designs and/or engine moding strategies result in 
unacceptable consequences as a result of modelling uncertainty. 
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V-D- Progress Toward AIPS Requirements Specification 

The HSV integrated trajectory and control analysis methodology developed under this task 
provides information that is essential to development of closed-loop G&C system design 
requirements. The demonstration winged-cone vehicle results have provided enough trajectory 
management strategy and performance sensitivity information (reviewed in the previous sub- 
section) to immediately support development of a strawman G&C requirements specification. A 
complete G&C specification would permit characterization of a closed-loop G&C algorithm, 
including input/output processing rates and redundancy requirements. Completion of all these 
steps, to accommodate assessment of expected overall G&C related throughput and compu- 
tational processing rates, is a mandatory part of the AIPS requirements specification [15]. 
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VI. Conclusions 


A generic methodology to permit integrated trajectory/control analysis, including the effects 
of physically derived constraints, has been successfully adapted for assessment of single- 
stage-to-orbit hypersonic vehicles with multi-mode propulsion systems. Both an overall trajectory 
evaluation algorithm and a preliminary near-fuel-optimal trajectoiy development tool have been 
mathematically formulated and implemented for use in a work station computational environment. 
The methodology is capable of assessing the relative performance of alternative vehicles and/or 
subsystem technologies with respect to a single mathematically valid measure. To assure good 
numerical behavior of the analysis algorithms, a cubic spline based data smoothing routine was 
developed and included in the resulting tool to permit use of multidimensional tabulated vehicle 
models. 

Hypersonic vehicle designs demand complex integration of the airframes and propulsion 
systems. As a result, there is an expectation that there will be a strong dependence between 
configuration details, the mission performance, and trajectory management and control strategies. 
To assure that an overall vehicle design would have a high probability of achieving mission 
objectives, it was thought to be essential to have a means to perform an early assessment of how 
specific vehicle designs and their guidance and control system characteristics would interact. 
Justification of this argument for the research demanded a specific demonstration of the generic 
analysis tool. 

In support of the demonstration objective, two applicable vehicle dynamics models were 
derived including all relevant analytic partial derivatives. A nine state model incorporates all 
translational and longitudinal (pitch) dynamics and uses throttle and elevon controls. A somewhat 
simpler seven state model should help bootstrap to the use of the nine state model. The simpler 
model restricts treatment of longitudinal effects to the calculation of elevon deflections and the 
resulting incremental forces which maintain vehicle pitch trim. Angle of attack in the seven state 
model replaces elevon deflection in the nine state model as a control, which should permit more 
intuitive analysis to be applied. 

A specific test case of the hypersonic vehicle integrated trajectory/control analysis meth- 
odology was defined, performed, and evaluated to demonstrate the power of the algorithm to 
derive information of critical importance to developers of guidance and control systems. The 
example vehicle model was based on an air-breathing accelerator, winged-cone vehicle database. 
Performance from horizontal takeoff to 20,000 ft/s relative airspeed on a polar orbit track was 
determined. A dynamic pressure constraint of 1000 psf was enforced. The following list 
summarizes some of the most significant conclusions: 

- Three near-fuel-optimal flight segments were derived that satisfy the applied constraints. 

They successively cover flight from horizontal takeoff to Mach 2, Mach 2 to 1 0,000 ft/s, 
and 10,000 ft/s to 20,000 ft/s (where velocity is specified in terms of relative air speed). 
Segment boundary condition matching criteria were enforced to permit the merger of 
the individual segment results into a single full trajectory usable for overall optimization. 

- The analysis methodology was able to identify previously unspecified discontinuities and 

propulsive phasing characteristics in the vehicle model. This information can be used 
to directly treat these strong model nonlinearities so as to reduce the computational 
burden in subsequent work. 

- Flight regimes with strong performance sensitivities and inadequate model fidelity were 

highlighted. Failure to consider these effects can leave the vehicle design at great risk 
of unacceptable performance losses. 

- Data to determine lower bound control effector response times and the remaining effector 

control authority after accounting for trim was generated. 
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- Drag losses due to required control effector activity was shown to be very high, demon- 
strating the importance of also including vehicle design characteristics in the integrated 
analysis. 

When analysis results are considered collectively, they provide much of the information 
needed to define guidance and control system design requirements. The methodology supports 
treatment of important effects not considered in the demonstration case, such as accommodating 
thermal constraints and determining efficient propulsive mode phasing, which would provide even 
more realistic inputs to the design process. 
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Appendix A 

Nine State Function Partial Derivatives for the Example Model 

A.1 . Derivatives with Respect to States 


The matrix f v for the example application with nine states (derived from differentiating Eq. 
V.19) is too large to print as a single equation, so it is specified here on a column by column 
basis. The equations for f v are simplified by the knowledge that the following derivatives have 


zero value: — 

d v 
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(Eq. A. 6) 


QCr Mjp. 



(Eq. A. 7) 


(Eq. A. 8) 



(Eq. A. 9) 



where: 


= sign\\)[~ sin t|i( F T cos y + l : N v,in y) — + cor>i|;l — — cosy* -j— siny 


+ cog ip — (-/•' 7 si n y + /•' N cos y ) ] 

b X 


dF 2 . by f t dF n . 

= sign cos r co;; y + /• N sin y) T“ * sm W T7 cos Y + -777 sm Y 

£ X (7 X \ (7 X O X 


+ simp — (- Fr sin y + / ; n cosy)] 
3 .v 


aF 3 bF T dF N dF c by c . 

— - -sin y + cosy + —— - — (F T co$y + F w sin y) 


b x b x 


b X b X 


>Q . ~ dC > 


- = cS reiyC + Q -^7 J “(Osina* /- cos a + <i> * 7 , si n 6 , ) 

(bD bL . ba\ 

-Xrr\ — s i n a + — cosa + (Dcosa- Z.sina) — 
CC \bx bx bxj 


— = f 0 0 0 


( *' 4 ) 2 + ( * 5 ) 2 ( * 4 ) 2 + ( * 5 ) ' 


0 0 0 0 


bF T bl bD da 

= — cos a — T sin a 

b X bx bx bx 

bF N bT . bL ba 

= sin a T cos a 

bx bx b X b X 


(Eq. A. 10) 

(Eq. A.1 1) 

(Eq. A. 12) 
(Eq. A. 13) 

(Eq. A. 14) 
(Eq. A. 15) 
(Eq. A. 16) 


52 


d I± 

dx 


/ 


dy^ 

dxj 


0 

0 

0 

-v 4 -v 6 


\ 


0'/) X n/(*4) 2 + (*5) : 


v\.v 


\ 


(i'/) 2 V( X 4) 2+ (-V 5 ) : 
~V(*4 ) 2+ C^s) 2 

O'/)’ 

0 
0 
0 


/ 


-2ii e .v, 

( x ,) 3 


0 0 0 0 0 0 0 




Ol)' 


3C 

L ' m 


m d C m 


(fc dC m<t 


cC 


m q d ( q 


d X 


dx 


d X 2 V f dx 


2 axV v 



dx 

da 

dx 

dM f 

dx 

»c m . 

dc m 

a d a 

dC <n 0 

0 

dMf 

dx 

da 

"dx + 

dM f 

dx 

d C m q 

dC.n, 

jd a 

dC 

dMf 

dx 

da 

dx 

dM, 

dx 



X 4 v 


■\ 8 


-.Vs.V 


S '- 8 


-X 6 X e 


(W/) 3 O,) 3 ("/)' 

= 100000000 


0 

dx 


- 0 


cc 


3x< 


(Eq. A. 17) 
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(Eq. A.26) 
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(Eq. A.35) 
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To complete the evaluation of f x , the following partial derivatives are derived from the 
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A.2. Derivatives with Respect to Controls 

All the equations required to compute the matrix / u for the example application with nine 
states, required for use in Eqs. V.49,50,54, are derived from differentiating Eq. V. 1 9. The equations 
for f u are simplified by the knowledge that the following partial derivatives have zero value: 
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where: 
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The matrix f x for the example application with seven states (derived from differentiating 
Eq. V.70) is too large to print as a single equation, so it is specified here on a column by column 
basis. The equations for f v are simplified by the knowledge that the following derivatives have 

zero value: i£ . . i£ . • L_i Simplifications are also made by utilizing the knowledge 
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(Eq. B.24) 


*Ci..dM f t *C L .db a 

dx d M f dx db e dx 

dC , dC T dM f dC T dQ 
~ ~ ~\ T f dx + ~dQ ~dx 


(Eq. B.25) 


dM f 

dx 


-v / da 
a A dr 


M, 

0 0 


a 
V\ 
av f 


l>2 
av f 



d6._pC„ t V '(( dC m% i dC m \dMj 

dx (_ db r J V V dM i + dM ! ) dx 


d f -V C q 

JZ\QcS raf 


( Dsin a + L cos a) 


<</ 


QcS, 


«/ 


dD d L 

— sin a + — cos a 

d x d x 


d f X cc ^ _ I f d x cc x CG d Q A 

^^QcS rp/ J"QcS re/ V~^ Q~dxJ 


dx CG 

dx 


0 0 0 0 0 0 


dx cc 
dm 


(Eq. B.26) 
(Eq. B.27) 


(Eq. B.28) 


(Eq. B.29) 


(Eq. B.30) 


To complete the evaluation of f x , the following partial derivatives are derived from the 
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. Also a and — are derived from a standard atmosphere model. 


B.2. Derivatives with Respect to Controls 

All the equations required to compute the matrix f u for the example application with seven 
states, required for use in Eqs. V.49,50,54, are derived from differentiating Eq. V.70. The middle 
column of / u only applies to leg 3 of the ascent trajectory: 
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(Eq. B.40) 
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Note that the middle terms in the above matrix equations for the partials with respect to 
controls disappear for legs 1 and 2 of the ascent trajectory. 

To complete the evaluation of f u , the following partial derivatives are derived from the 
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* + a i4, a dU 3 3u 3 iu 3 iu 3 
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Xcg(ft) 


1.4 1.6 1.8 2 2.2 2.4 2.6 2.8 3 

mass in lbs xl05 

Figure V-2. Vehicle Mass Model Curve Fit 

Xcg fit 



mass in lbs xl05 


Figure V-3. Vehicle Pitch Inertia Moment Model Curve Fit 
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X - ' Y - Z : Earth centered inertial frame 

x - y - z : Earth fixed frame (rotates with Earth) 

1 - 2 - 3 : Vehicle fixed LVLH frame 


Figure V-4. Coordinate Systems for Dynamics Equations 




67 


Figure V-5. Integrated Trajectory/Control Optimization Algorithm Structure 
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Figure V-9. Demonstration Case - Equivalence Ratio vs Time 
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Figure V-10. Demonstration Case - Elevon Deflection vs Mach Number 
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Rgure V-1 1 . Demonstration Case - Elevon Deflection vs Time 
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Figure V-13. Demonstration Case - Dynamic Pressure Derivative vs Mach Number 
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Figure V-14: Demonstration Case - Flight Path Angle vs Mach Number 
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Rgure V-17. Demonstration Case -Inertial Acceleration vs Mach Number 
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Figure V-19. Demonstration Case - Mass Flow vs Mach Number 
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